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AERODYNAMIC DESIGN AND ANALYSIS SYSTEM 
FOR SUPERSONIC AIRCRAFT 

PART 1 - GENERAL DESCRIPTION AND THEORETICAL DEVELOPMENT 

W. D. Middleton and J. L. Lundry 
Boeing Commercial Airplane Company 

. 1.0 SUMMARY 


An integrated system of computer programs has been developed for 
the design and analysis of supersonic configurations. The system 
uses linearized theory methods for the calculation of surface 
pressures and supersonic area rule concepts in combination with 
linearized theory for calculation of aerodynamic force 
coefficients. Interactive graphics are optional at the user's 
request . 

The integrated system consists of an executive "driver” and seven 
basic computer programs including a geometry input module, which 
are used to build up the force coefficients of a selected 
configuration. 

The description of the design and analysis system is broken into 
three parts: 

Part 1 - General Description and Theoretical Development 

Part 2 — User's Manual 

Part 3 - Computer Program Description 

This part presents a general description of the system and 
describes the theoretical methods used. 



2.0 INTRODUCTION 


Over a period of years, NASA-Langley has developed a basic 
computerized series of supersonic design and analysis methods for 
aerodynamic configuration studies (reference 1). The methods are 
characterized by their reliability in use and input simplicity. 

The Boeing Company has extended this basic series of methods and 
combined them into an integrated system of computer programs. The 
extensions to the methods provide several new features: 

• Addition of a near-field (thickness pressure) wave 
drag program, to complement the existing supersonic 
area rule program. 

• Improved modeling of fuselage in lifting surface 
design and analysis programs. 

• Addition of configuration-dependent loadings in wing 
design program, so that the wing design is performed 
in the presence of fuselage and nacelle effects. 

• Addition of pressure limiting terms in the lifting 
pressure programs, to constrain the linear theory 
solution. 

• Optional CRT displays of selected program input and 
output data and provisions for limited user editing 
and intervention. 

A plot module is included in the system to produce configuration 
drawings, and a common geometry module is used to permit a single 
geometry input for all programs. 

The basis of the system is supersonic linearized theory, modified 
in two respects: 

• The "Whitham" correction to disturbance positioning 
is used in the propagation of body pressure fields. 

• The wing lifting pressure modules contain an 
optional limiting pressure feature to control the 
permissable level of upper surface pressure 
coefficient. 


Superposition is used to build up the theoretical force 
coefficients of a selected configuration. 

The goals of the integrated system have been to develop an easily 
used supersonic design and analysis capability, with recognition 
of the need for constraints on linear theory methods to provide 
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physical realism, and with inclusion of interactive display for 
increased design control over optimization cycles. 
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3.0 SYMBOLS 


A 

A. 

1 

A (x) 

B 

c 

C 

c 

C D 

C F 

C L 

C 

m 

C 

no 

C 

P 

AC 

P 

D 

F(y) 

h(Z) 

k 

L 

M 

N 

P 

g 

R 

B 


element fraction 
load strength factor 
cross sectional area at x 
trailing edge element fraction 
pitching moment reference length 
wing tip element fraction 
chord 

drag coefficient, D/qS 
skin friction coefficient 
lift coefficient, lift/qS 

pitching moment coefficient, pitching moment/qSc 

pitching moment coefficient at zero lift 

local pressure coefficient, (p - p M )/q 

lifting pressure coefficient (lower surface pressure 
coefficient minus upper surface pressure coefficient) 

drag force 

Whitham function 

decay function for F(y) calculation 
Mach number parameter 

grid element variable in x ( lengthwise) direction 
Mach number 

grid element variable in spanwise direction 
pressure 

dynamic pressure, .7 p M m 2 
body radius 
influence function 


r 

s 

T 

t 

0 

u 

W(Z) 

w 

X 

y 

z 

0 

a 

e 

B 

5 

A 

n 

P 

♦ 

0 

T 

M 

Sub sc ri 


radial distance to point in flow field 
reference area or cross-sectional area 
tern perature 

variable of integration (figure 4.3-4) 
free stream velocity 
x perturbation velocity 
body pressure decay function 
upwash 

lengthwise variable 

span station or distance variable 

vertical variable 

decay function parameter 

angle of attack 

local fuselage camber angle 

Mach number parameter = yM 2 - 1 

x variable of integration 

surface shape (equation 13) or Lagrange 
(figure 4.4-4) 

spanwise variable of integration 
density 

velocity potential 
radial angle variable 
region of influence 
dynamic viscosity coefficient 
pts: 


(figure 4.3-3) 


multiplier 


00 


6 


i 


free stream conditions 
incompressible condition (equation 5) 


i, j loading numbers (equation 31) 

t tilted lengthwise variable (figure 4. 3-5) 

Superscripts: 

i reference temperature condition (equation 2) or 

first derivative ’with respect to x 

.. second derivative with respect to x 

* field point grid element 
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4.0 DISCUSSION 


The integrated supersonic design and analysis system is shown 
schematically in figure 4.0-1. It consists of an executive 
"driver" and seven basic computer programs. The individual 
programs, or modules, provide data for configuration design or 
analysis as follows: 

• Skin friction is computed using turbulent flat plate 
theory . 

• Have drag is calculated from either far-field 
(supersonic area rule) or near-field (surface 
pressure integration) methods. The far-field 
method is used for wave drag coefficient calculations 
and for fuselage optimization according to area 

rule concepts. The near-field method is used 
primarily as an analysis tool, where detailed 
pressure distributions are of interest. 

• Lifting pressures, drag-due-to-lif t, pitching moment 
and trim drag are computed from the lift analysis 
program, which breaks arbitrary wing/fuselage/ 
canard/nacelles/horizontal tail configurations into a 
mosaic of "Hach-box" rectilinear elements which are 
employed in linear theory solutions. A complementary 
wing design and optimization program solves for the 
wing shape required to support an optimized pressure 
distribution at a specified flight condition. 

• The plot module draws configuration pictures according 
to size and view inputs. 

• The geometry module handles configuration geometry 
for the system. 

« Interactive graphics, for data display or editing, 
are optional in the system at user's request. 

Operating in the analysis mode, the force coefficients of a 
selected configuration are built up through superposition as shown 
in figure 4.0-2. 

The entire design and analysis system is a single overlaid 
program. The executive level of the system controls module 
execution by means of special identification cards in the input 
data. Transfer of data in the system between modules is handled 
by disk storage and common blocks. 

All configuration geometry is read (or updated) in the geometry 
module. Geometry inputs are patterned after those of the NASA-LRC 
configuration plotting program. All "paneling" of the 
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SUPERPOSITION METHOD OF DRAG ANALYSIS 



ZERO LIFT WAVE DRAG 


FIGURE 4. 0-2. -DRAG BU/L DUP 
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configuration for theoretical analyses is accomplished within the 
programs, and the user prepares only "drawing-type" geometry. 

Additional data on system structure and input and output formats 
are presented in parts 2 and 3. The theoretical methods used in 
the modules are described in the remainder of this document. 


4.1 Skin Friction Drag 


The skin friction drag of an arbitrary configuration is calculated 
from the t' method of reference 2. Smooth flat plate, adiabatic- 
wall, and turbulent boundary layer conditions are assumed, with 
transition assumed to occur at the leading edge of each 
configuration component. 

The theory and experimental verification of the T’ method are 
given in detail in reference 2. The essential elements of the 
method are presented in Appendix C of the reference, and* are 
summarized here for completeness. 

The T' method is based on the calculation of a compressible skin 
friction coefficient, C p , from a reference skin friction 
coefficient, C p , for a selected Mach number, M,*,, Reynolds number, 
Roo, and adiabatic wall temperature, T w . (The subscript «> is used 
throughout to denote free stream conditions) . 

The wall temperature ratio, T^j/T^ , is calculated from one- 
dimensional relationships assuming a wall recovery factor of .89. 


= 1 + 0.178 Mi (1) 

1 OO 

The prime superscript (T’ , R' , etc.,) refers to conditions at 
which incompressible flow relations must be evaluated in order to 

in reference 2 


( 2 ) 


The Reynolds number relationship is 

A L _ i ' (3) 

R ~ (fj (t) 

where the viscosity ratio is given by the Sutherland equation 


represent compressible flow/ 
obtained the key relationship 


Sommer and Short 


T ' 
T_ 


= 1 + 0.035 M „ + 
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y' / T oo + 216 y 

Ko ~ TO NT' + 216 / 


with T in degrees Sankine. 

The incompressible skin friction coefficient is given by the 
Karman-Schoenherr equation 


which gives the 
analogy: 


This equation is solved iteratively for Cp' . The desired 
compressible skin friction coefficient is obtained from C ' by 

c 

°F = V ( T oo /T ' ) (7) 

Most of the skin friction program is involved with computing 
wetted areas and reference lengths for the configuration. 
Components which may exhibit large variations in reference length 
(such as the wing, tail, canard etc.,) are broken into strips in 
order to calculate Cp accurately. 


0.242 

10610 


( c ?i »i) 


(5) 


corresponding relationship based on the 


( 6 ) 


4.2 Far-Field Have Drag 


The far-field wave drag program computes the zero lift wave drag 
of an arbitrary configuration by means of the supersonic area 
rule. A description of the program has been given in reference 3, 
and is summarized here. The method assumes not only an integrated 
total lift of zero, but local lift everywhere zero, and thus 
includes none of the wave drag contribution associated with the 
generation of lift. 

The supersonic area rule is a generalization of the transonic area 
rule, which states that the transonic wave drag of an aircraft is 
essentially the same as the wave drag of an equivalent body of 
revolution having the same cross-sectional area distribution. In 
the supersonic area rule procedure, several equivalent bodies of 
revolution are produced by passing a series of parallel cutting 
planes through the configuration, as shown below. The cutting 
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planes are inclined with respect to the aircraft axis at the Mach 
anqle M , and a single equivalent body is produced for the series 
of cutting planes at a constant azimuthal angle, 0 . The area of 
the equivalent body at each station is the projection of the area 
intercepted by the cutting plane onto a plane normal to the 
aircraft axis. 



The wave drag of each equivalent body is calculated from the von 
Karman slender-body formula (reference 4) , which gives the drag as 
a function of the free stream conditions and equivalent-body area 
distribution. 


pu 2 r ^ 

D(6) = I I A " A " ^ 106 I X 1 " X 2^ **1 ^2 

where Xi and X 2 are lengthwise variables of 
integration and A" is the second derivative of the 
body area distribution. 


The wave drag of the aircraft at a given Mach number is calculated 
from the integrated average of the equivalent-body wave drags. 


D = 


_1 

2ir 



2tt 

D(0) d6 


(9) 


A useful feature of the supersonic area rule occurs in the 
optimization of a fuselage area distribution to minimize the wave 
drag of a wing-fuselage combination. According to Sheppard 
(reference 5) , the wave drag of a wing-fuselage can be written as: 



wing 


(t) -(i) 

wing eq. body wing eq. 

+ body body 


( 10 ) 
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where the wing equivalent body is a body of revolution obtained 
by averaging from 0 = 0 to 2n all the projected areas intercepted 
by cutting planes passing through the wing for each X station. 

Considering a wing- fuselage , then, and assuming that the wing 
geometry is fixed, both the first and third terms of eguation 10 
are fixed. This leaves optimization of the second term only, so 
that to design the minimum drag configuration, the fuselage must 
be contoured to produce a wing-eguivalent-body-plus-f uselage 
having the shape of a minimum drag body of revolution. 

The wave drag program uses the Eminton-Lord fairing through a 
discrete set of points, which defines the shape of a minimum drag 
body through the points. Using this fairing, and identifying a 
few "control-points” on the fuselage area definition, an optimized 
fuselage area definition may be obtained as the difference between 
the wing equivalent body and the combined wing-equivalent -body- 
plus- fuselage. 

The wave drag program has mechanized this solution, and also 
includes the effects of nacelles, tail, etc., which are included 
similarly. The shape of the minimum drag fuselage area 
distribution, and the drag of the complete configuration including 
the minimum drag fuselage, may be calculated for each 
configuration input. 

Because of the configuration generality that can be handled, the 
far-field program is the primary source of zero-lift wave drag 
coefficients in the design and analysis system. A complementary 
near-field program, more restricted in scope but which computes 
detailed thickness pressures, is described in the following 
section. 


4.3 Near-Field Wave Drag 


Isolated Component Pressures 


The near-field wave drag program computes zero-lift thickness 
pressure distributions for an arbitrary wing -body-nacelle 
configuration in supersonic flow. The pressure distributions are 
integrated over the cross-sectional areas of the configuration to 
obtain the resultant drag force. Three basic calculations are 
performed to obtain the required pressure fields: 

• Thickness pressure distribution for a wing of 
arbitrary section and planform. 

• Thickness pressure distribution on surface of . 
fuselage or nacelles. 
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'•Shitham'* near-field calculations to define pres- 
sure distributions propagating from fuselage or 
nacelles. 


Superposition is used to calculate the interference drag terms 
associated with the pressure field of a given component acting on 
the surfaces of the other components. 


The surface pressure distributions calculated include some (from 
the nacelles) which contribute lift due to a non-symmetrical 
distribution of volume, so that the drag calculations which 
include these loadings contain wave drag due to lift 
contributions. However, by selection of program options and 
subtraction of certain drag contributions, a zero lift wave drag 
comparable to the far-field wave drag value may be obtained, as 
described on page 37. 


Wing thickness p re ssures. - The surface pressure coefficient on 

the upper (or lower) surface of a flat-mean-line wing of 
symmetrical surface shape is obtained by first calculating the 
corresponding velocity potential, differentiating with respect to 
x (to get u) , and then computing the pressure coefficient from the 
linear theory approximation, C p = -2U/UO,,. 


The velocity potential computation, from reference 13, is: 


where 


♦ (x,y) = - — 




X dn d£ 


^(x-C) 2 - g 2 (y-n) 2 


<t> (x, y) 
X 


0 

t 

V 


velocity potential at a defined wing field 
point (x,y) 

surface slope (dz/dx) of wing section at 
a wing i ntegration point 
VM 2 - _ 1 

x variable of integration 
y variable of integration 

subscript denoting interval of integration 
(surface of the wing planform within the 
Mach forecone from x,y) 


( 11 ) 


The wing thickness pressure coefficient is: 

. e p (x, y ) = ^=- 2 ^=-2ii <«) 

where ? = local pressure at x,y 

= free stream static pressure 
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u = x perturbation velocity 
Uoo»q<» = free stream velocity and dynamic pressure 

Integration of the velocity potential equation is performed by 
representing the wing as a grid of rectangular elements, and 
substituting a numerical summation for the integration. The wing 
element "fiach-box M system exactly follows the arrangement 
described in Section 4.4 (for the lifting pressure case) and is 
summarized in figure 4.3-1. The grid elements, identified by L 
and N, are defined such that L is numerically equal to x and N is 
numerically equal to , where x and p y take on only integer 
values. Partial grid elements along the wing leading and trailing 
edges are used to permit a closer approximation to the actual wing 
planform. The grid system of figure 4.3-1 is rather coarse for 
illustrative purposes; in actual usage, many more grid elements 
are employed. 

With respect to a specified field point x,y, the upstream region 
of influence t is approximated by the shaded grid elements in 
figure 4.3-1. Each of the shaded elements has associated with it 
an influence factor, ¥, which relates the effect of the element 
and its surface slope to the total velocity potential at the field 
point. 


The factor B is obtained by an approximate solution to Equation 11 
as follows: 


Let y- n = u so that du = -dn 


A / x 1/’/’ Adudt (13) 

’'"’’■’W 77? 

Substituting an average value of (x -£) = (x -£) and integrating 
with respect to u. 


* U ’ y) *7 [ 

Restoring y - n.= u 

r T 

, / x 1 I , . -1 B(y-n) 

* (x - y) = 7TL Sln trtJ 

n l 

Eor a given grid element in the Mach forecone, 
interval to n 2 gives 

0(y~n. ) = N* - N ♦ 0.5 

0(y-n 2 ) = N* “ N -o.5 


(HO 


(15) 


the integration 
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L* 


L + 0.5 


Also, X-£ 

L* and N* define the field point grid element and L and N define 
an element in the forecone. The numerical summation then becomes 


<t> (x,y) 


= X L,N 

T 



N* - N - 0.5 
L* - L + 0.5 


sm 


-1 H* - K 


L* - L 


± 0.5 1 
+ 0.5 J 


(16) 


with R being the quantity between the brackets for the element 

(L,N). 


The character of the R function is such that the sum of R for all 
elements across the Mach forecone at a constant L value is equal 
to -ir . in considering the case of an infinite rectangular wing 
of constant thickness shape (wedge section) , it is seen that the 
-if summation must be the case. Since the <t> function for 
successive L* rows is increased by 

L=const . 


so that 


£-5- 


A 

P 


and 


C = 
P 


-2 


9 X 
to " 2_ b 


which is the established linear theory result, 
the Mach forecone from a typical field element 
figure 4.3-2. 


Values of R within 
are presented in 


The computer program for the wing thickness pressure calculations 
sets up the wing grid system, defines the surface slope X for 
each grid element, performs the summation for <f> at selected field 
elements, and differentiates 4> numerically to obtain C p . Have 
drag is obtained by summing the pressure coefficients times the 
corresponding surface slopes of the wing. 


In the wing surface slope definition, slopes of the elements 
occupying space covered by the fuselage (if there is one) are 
set equal to zero to eliminate their drag contribution. Partial 
elements are used at the wing-body intersection to improve the 
geometry definition. 


As was found to be the case in lifting surface programs involving 
the grid system used in the wing calculations, some oscillation in 
computed pressure coefficient values occurs. This has been 
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suppressed by a smoothing of the calculated pressure coefficient 
values, using a 5 term equation of the form; 


C = 
P 


• 333 v-2 y, + - 66T v-i y. 


L*-l ^ 

.333 a Ij *_2 + \,*-l + 


+ 1.0 


(17) 


where A = factor defining size of associated element (1 if whole 
element, fraction if not) 


Fjis.ela_q^e and nacelle thickness pressures . - The pressure 

distribution on the surface of fuselage or nacelles is obtained 
from a method developed at the Boeing Company by R. M . Kulfan, 
based on the Lighthill theory of reference 6. The method is 
applicable to bodies having either smooth area distributions or 
bodies with slope discontinuities. As used in the near-field 
program, smooth area distributions are assumed except (if 
required) at the nose or aft end of the body. Open-nose bodies, 
such as nacelles, are permissable. 


The solution technique is summarized in figure 4.3-3. An axial 
perturbation velocity is calculated which is a function of the 
body cross-sectional area growth (and radius, distribution) and a 
decay function which relates area growth to its effect on a given 
field point. 


The axial perturbation velocity, u, is given by 


u(x) = 


- 4 . 


W(z) 


e r 


d S' 


(K) 


(€) 


(18) 


where 


x 

W(Z) 

£ 

z 

(£) 


body field station 
decay function 
x variable of integration 
position function = x ~ £ 

body radius at £ vw 

first derivative of body cross-sectional 
area (S) at £ 


The pressure coefficient is calculated from one-dimensional flow 
relationships (rather than the simpler linear theory 
approximation), as shown in figure 4.3-3. The decay function, 
W (Z) , is tabulated on page 23. 


fuse la ge and nacelle flow-field pressures. - The pressure field 
propagated by the fuselage or nacelles is calculated using the 
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DECAY FUNCTION 
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FIGURE 4.3-3.— BOD Y SURFACE PRESSURES (LIGHTHILL THEOR Y) 


DECAY FUNCTION 
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Rhitham solution of reference 7. This solution, which is 
converted to a method amenable to the digital computer in 
reference 8, calculates shock wave positions and strengths through 
a modification of linear theory results. 


The Whitham solution begins with the calculation of a function, 
F(y), which is dependent on body geometry: 


where y 

h (Z) 

t 

z 

s'(t) 

^end 


F(y) 


■*/ 


^end 


| 

2 

J 

0 

htjj 

station 


_1_ 

2 


h ( Z ) d S' 


(t) 


(similar to body pressure 


body fie] 
decay function 
calculation) 

y variable of integration 
t position function 

first derivative of body cross-sectional 
area (S) at t 

y at end of integration (see figure 4.3-4) 


(19) 


The solution for F(y) is very much like that of the body thickness 
pressure solution, with the exception that the integration is 
carried out to the point (y erM j ) at which the aft-running Mach line 
from the body centerline at y emerges from the body contour. The 
F(y) calculation is summarized in figure 4.3-4, and the decay 
function h (Z) is tabulated on page 27. 


Development of the near-field pressure signature at a radial 
distance (r) from the body centerline is summarized in figure 4.3- 
5. The F(y) function is tilted, as shown, by displacing the F(y) 
function according to 


where 


Ay = k/F F(y) 


( 20 ) 


. _ 1.6973 M 
' 6 1 - 5 


This process results in a double-valued pressure signature where 
shock waves are present, . a situation which is corrected by an 
area -bala ncing technique. The area-balancing (shaded lobes are 
equal in area) defines shock wave locations and strengths. The 
resultant tilted/area-balanced signature is then converted to 
pressure coefficient by 


= p ~ p « „ 2 F(y) ( 21 ) 

P 3co V 2f5r 
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The position of the pressure signature at radius (r) is then given 
by 

x = y t + Br (22) 

where x is the longitudinal distance from the body nose. The 
"tilting" thus produces a correction in pressure signature 
positioning, relative to simple translation along a Mach line, 
which is the result of the remarkable theory of reference 7. 


Interference Calculations 


The previous sections summarized the basic pressure calculations 
performed by the near-field wave drag program. Integration of the 
wing and fuselage or nacelle pressure distributions over their 
respective surfaces gives the "isolated" wave drags of the 
components. Interference drags are obtained by superposition; 
i.e., the pressure field of each component is imposed upon the 
surfaces of the other components to calculate the resultant 
interference forces. 


The computer solution allows for up to 3 pairs of nacelles located 
external to the wing-fuselage (or 2 pairs plus a single nacelle at 
Y=0) . The nacelles may be either above or below the wing (or 
both) . 


Wing -fus elage interfe ren ce. - Fuselage-on-wing interference is 

obtained by calculating the near-field pressure signature from the 
fuselage at selected spanwise stations, and imposing them upon the 
corresponding wing sections as a bouyancy field. The spanwise 
stations are the same as those used in the isolated wing thickness 
pressure calculations. 


Wing-on-fuselage interference is obtained by computing wing 
thickness pressures in the area occupied by the fuselage, after 
the wing surface slopes in the fuselage area have been set egual 
to zero. The thickness pressures thus calculated are "carry-over" 
pressures, and are imposed upon the fuselage surface slopes by 
transposing them aft along Hach lines to the fuselage. The 
fuselage area covered by the wing root is deleted from the 
interference term. 


Typical results from wing-body calculations are presented in 
figure 4.3-6. The wing data, both isolated and f uselage-on-wing 
interference, are converted to section drag form. 
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FIGURE 4.3-6.- WING-BODY SOLUTION, M = 2.6 







Nacelle interference terms are calculated like the wing- body 
interference terms, with the provision that two alternate 
solutions may be specified for the nacelle-on-w ing interference 
term. Available experimental data do not make it clear whether a 
"wrap" or "glance” solution, as shown in figure 4.3-7, is more 
correct. Since the nacelle- on- wing interference term is 
substantial, both solutions are available in the program 
(controlled by an input code) . 

The nacelle interference terms are calculated as follows: 

I§£ell§zon-wing_term. - The nacelle pressure field acting on the 
wing is obtained by calculating nacelle pressure signatures at the 
same spanwise stations used for the wing thickness pressures (plus 
extra stations immediately adjacent to the nacelle centerlines) , 
then defining a composite signature by summing together the 
contributions from all nacelles. The nacelle pressure 
coefficients are doubled to account for reflection from the wing 
surface. 

In the case of nacelles both above and below the wing, separate 
nacelle pressure fields for nacelles below, and then above, the 
wing are calculated and then summed. 

A unique feature of the near-field approach lies in the solution 
choices available, since the pressure signatures generated by a 
nacelle to act on the wing surface may be terminated on 
encountering another nacelle or may' be allowed to pass 
undiminished around (or through) it. In the first case, the 
"glance" solution, the nacelle pressure field is assumed to 
reflect from other nacelles in such a way that the reflected field 
exerts no further influence upon the wing surface. In the second 
case, the "wrap" solution, the nacelle pressure field is assumed 
to propagate around the other nacelles as if no obstruction were 
offered, i.e., as if they were transparent to the pressure field. 

A comparison of the theoretical nacelle pressure field with 
experimental data is presented in figures 4.3-8 and 4.3-9 at Mach 
numbers of 2.7 and 1.1. At M =2.7, there is no difference 
between wrap and glance solutions. At M = 1.1# the solution 
choice appears to favor the glance approach except on the outboard 
wing. However, the experimental pressure points are too sparse to 
permit firm conclusions. 

Nacell e- on- fu se l age term . - The effect of the nacelles on the 

fuselage area distribution is built up by integrating each 
nacelle,* s pressure signature upon the fuselage area growth, with 
the provision that nacelles below the wing affect only the 
fuselage area below the wing in the region of the wing (and vice- 
versa for above-wing nacelles). In the fuselage region below the 
wing, or above the wing, the nacelle pressures are doubled to 
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FIGURE 4.3-7.-NACELLE GENERA TED WING PRESSURE FIELD CONCEPTS 









account for reflection (as is done in the nacelle pressure field 
case) . 


Fuselage -on -nacelle term. - The fuselage pressure signature is 

imposed on each nacelle surface as a bouyancy force. 


Wing-on- na cel le_term. - Wing thickness pressures at span stations 
where nacelles are located are transferred aft along Mach lines 
from the wing surface to the nacelle centerline to obtain the 
buoyancy field acting on each nacelle. 


Nacelle^on-nacelle_t er m. - The interference term of other nacelles 
acting on a selected nacelle is calculated by building up the 
composite buoyancy field, and then imposing it upon the nacelle 
surface. In the case of nacelles on opposite sides of the wing, 
the pressure signatures are cut off where intercepted by the wing. 

The nacelle-generated pressure field is assumed to pass 
undiminished around or through other nacelles that may be present 
when a particular nacelle pair interference term is being 
calculated much as in the "wrap" solution employed for nacelle on 
wing terms. The "wrap" or "glance" option is not provided in this 
case because of the generally negligible difference in results. 

Nacelle im ag e ef fe cts. - If the nacelle is located next to the 
wing (or body), an "image" nacelle is used to create the reflected 
buoyancy field. The reflected field is cut off forward or aft of 
the wing if the nacelle is not entirely under the wing. 

The principal image effect is caused by the nacelle pressure field 
reflecting off the wing back onto itself. However, the reflected 
nacelle pressure field acting on other nacelles is also computed. 
This solution also utilizes an image nacelle representation, and 
the reflected signature is restricted to whatever part of the 
"real" signature encountered a reflecting surface. 


Composite Thickness Pressure Signatures 
and configuration Drag 

Typical theoretical thickness pressure distributions acting on the 
configuration components are presented in figure 4.3-10; The net 
effect of these pressure signatures on drag is obtained by summing 
the various inputs and integrating them over the surface of the 
configuration. The corresponding drag terms are summarized in the 
data of figure 4.3-11. 
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FIGURE 4.3-10.-NEAR-FIELD PRESSURE SIGNATURES, M = 1.1 "WRAP” SOLUTION 
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FIGURE 4.3-11. -TYPICAL WA VE DRAG COEFFICIENT SUMMARY 
NEAR-FIELD PROGRAM, M= 1.1 


The near-field wave drag program has three calculation features 
that have no direct counterpart in the far-field wave drag 
program. These are: 

1) Wrap or glance solution for nacelle pressure field acting on 
wing. The wrap solution is nominally the same as the far- 
field solution. 

2) Nacelle image effects. Because of the transparency 
assumption of the far-field program pressure propagation, 
nacelle pressures do not reflect off adjacent surfaces. 
(Through addition and subtraction of separate calculations 
involving "image" nacelles, however, comparable drag data can 
be generated using the far-field program.) 

3) Above- or below- wing fuselage area separation in the 
n acelle-on- fuselage term. Directly comparable results are 
obtained only for mid-wing configurations. 

Using the "wrap" solution, a mid-wing arrangement, and subtracting 
the nacelle image drag terms, a direct comparison between near- 
fiald and far-field program drag calculations can be made. 


4.4 Drag-Due-to-Lift (Design and Analysis) 

The wing design and lift analysis -modules are separate lifting 
surface methods which solve the direct or inverse problems of: 

• Design - to define the wing camber surface shape 
required to produce a selected lifting pressure 
distribution. The wing design program includes methods 
for defining an optimum (least drag) pressure 
distribution. 

• Lift analysis - to define the lifting pressure 
distribution acting on a given wing camber surface, and 
calculate the associated force coefficients. 


Originally, these programs were developed to utilize the "Mach- 
box" method of wing representation and evaluation of linear theory 
integral equations (references 9 through 11). In the Mach - box 
method, the wing is replaced by a grid of small rectilinear 
elements. Since many elements (in the thousands) can be used, a 
detailed description of complicated surface shapes, with 
associated computational accuracy, is possible. 


The wing design and lift analysis programs have been expanded 
several times to add features (e.g., reference 11) and were 
recently reviewed in reference 12. The discussion of the 
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aerodynamic theory of the Mach-box method in reference 12 is quite 
detailed, and is summarized here for completeness. The numerical 
method for the "direct" case of wing design and optimization is 
given first, followed by the inverse case of lift analysis on page 

50. 


Design Case 


Camber surface for a given loading. - A typical wing planform 
described by a rectangular cartesian coordinate system is 
illustrated in figure 4.4-1. The y co ordinat e has been multiplied 
by the Mach number parameter, g = y ti 2 - i for convenience in 
mathematical manipulations. 

For a wing of zero thickness lying essentially in the z=0 plane, 
linearized theory for supersonic flow defines the wing surface 
shape necessary to support a specified lift distribution by the 
inteqral equation 

(23) 


az 

5^U,y) 


-A 


4 V x - y) + TS 



(x - £) ac (s,n) dn 
E 

(y - n) 2 V (x - £) 2 - 8 2 (y - n) 2 


which is slightly modified form of equation (77a) of reference 13. 
With respect to a specified field point x, y, the upstream region 
of influence, r , enclosed by the Mach forecone is shown by the 
shaded area. 


Equation 23 may be rewritten into an influence function form by 
introducing the factor R, such that 

^ (x,y) = =jj-AC p (x,y) j R(x-e,y-ri) AC p U,h) <*Bn (2U) 

where I is defined as 


H(x-5 ,y-n) 


x — e 

g 2 (y - n) 2 *\/(x - O 2 - e 2 (y - n) 2 


(25) 


- 8 is a function relating the local loading at point £ ,n to its 
influence in determining the downwash (or upwash) at the field 
point. 
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FIGURE 4.4 - 1. -CARTES! AN COORD! NA TE SYSTEM 


For the numerical evaluation of this integral, the wing is 
replaced by the grid system of figure 4.4-2. The grid elements, 
identified by L and N, are defined such that L is equal to X and 
N is equal to 0y , where X and /3y take on only integer values. 
Partial grid elements along the wing leading and trailing edges 
are used to permit a closer approximation to the actual wing 
planform. The grid system of figure 4.4-2 is rather coarse for 
illustrative purposes, in actual usage, many grid elements are 
employed. 


In the grid system, the field point x, y is located at the aft 
midpoint of a field point element L*, N*, and the region of 
integration, r , is approximated by the shaded grid elements. 


The integration required by equation 24 is performed for each 
element within the Mach forecone, considering A,c p constant over 
the element and using an average value of X- £ = L*-L*5. The 
resulting expression, derived in reference 12, is the factor S 
which is the value of the influence function for an element L, N: 


n/' T « r Trt V(L* - L + 0.5) 2 - (N* - N - 0.5) 2 
R(L*-L,N*-N) - (L * _ L +0 /5)(N* - N - 0.5) 


- V 


(L* - L + 0.5) 2 - (N» - N + 0.5) 2 
(L* - L + 0.5)(N* - N + 0.5) 


(26) 


The integral equation 23 may then be replaced by the summation 
given below, where the summation includes all elements within the 
forecone and on the wing planform. The factors A, B, and c are 
element fractions for the wing leading edge, trailing edge, and 
wing tip, respectively, to allow for partial elements at those 
locations. 


t- S -(L»,H») -4-AC (L»,N») 

* — _ P (27) 

6 >.R(L*-L,N*-N) A(L,N) B(L,N) C(L,N) AC (L,N) 

TH P 

T 


‘The character of the R function is such that, for_a given L*~L set 
of elements within the forecone, the sum of the fi values is zero, 
the single negative value at N*-N=0 balancing all the others (see 
figure 4. 4- 3) . _At L* =L, where there is only one element in the 
summation, the R value is zero. 


The physical significance of this R variation is that (for 
positive lift), all elements directly ahead of the field point 
element contribute only downwash and all other elements contribute 
upwash. An element at the leading edge near the wing tip of a 
subsonic leading edge wing, therefore, sees a concentrated upwash 
field. It is this upwash field that makes the subsonic leading 
edge twisted and cambered wing attractive from the standpoint of 
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FIGURE 4.4-3. — tNFL UENCE FACTOR CHARACTERISTICS 


drag-due-to-lift , since a local element may be inclined forward to 
produce both lift and thrust. Conversely, an element located at 
the trailing edge and centerline of a diamond planform produces 
lift inefficiently, since the element sees an influencing region 
that contributes downwash. 


In computing the local surface slope required to obtain a 
specified lifting pressure coefficient, the ACp value is known 
everywhere on the wing from loading formulas. In the wing- design 
computer program, the wing camber surface shape is found by 
integrating the surface slopes for all chordwise elements at 
selected span stations. Section values of drag (ACp times slope) , 
lift and pitching moment are then integrated spanwise to obtain 
wing force coefficients. 


O p timum. combination o f loadings. - Equation 27, which defines the 
surface slope to support a specified lifting distribution, 
requires a description of the design pressure distribution. This 
is obtained from Lagrange's method of undetermined multipliers, 
which provides a means of selecting an optimum combination of 
component loadings, yielding minimum drag subject to various 
imposed constraints on the allowable pressure distribution.. 


A" number of design point options are provided in the wing design 
program, involving the selection of loadings to be combined and 
the constraints desired. (The options are controlled by input 
codes, as described in the user's manual, part 2). 

The usual design conditions specified are: 

• Drag-due-to-lift of the wing be minimized at a 
given design lift coefficient, subject to an 
optional pitching moment constraint. 

• Effects of fuselage and nacelles be included in the 
design solution. 

• Constraints be applied to the design pressure 
distribution to provide physical realism. 

The design point options are discussed on page 49. 

The optimum loading selection is an extension of the methods used 
in references 9 and 11, through the addition of the configuration- 
dependent loadings (due to fuselage and nacelles) and the addition 
of the pressure constraint formulation. Also, more component 
loadings were added to provide flexibility in rearranging pressure 
distributions at the design point in order to meet the pressure 
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constraint condition. The maximum number of loadings is 17, 
defined in the table on page 45. 


The wing design method actually consists of four solution steps: 

• Flat wing solution 

• Calculation of aerodynamic characteristics of 

selected component loadings. 

• Optimization of loading combination 

• Camber shape calculation 


Flat_w ing_so lution . - The first step in the wing design process is 
to obtain the flat wing solution for the given wing planform and 
Mach number, using the analysis form of the lifting pressure 
calculation described on page 50. This solution is obtained to 
locate the wing aerodynamic center for use in computing the 
pitching moment at zero lift for the individual loadings; and 
also, because its drag-due-to-lift factor is an often-used 
reference point in judging twisted and cambered wing designs. 


Com po ne nt loadings c haracteri st ics. - For each selected component 
loading, section lift, drag, and pitching moment coefficients are 
evaluated at a series of span stations and integrated spanwise to 
obtain overall lift, drag, and pitching moment. In addition, the 
interference drag coefficients between loadings are computed, 
where the drag interference between loadings i and j is given by 

»N* ) , 

(28) 

+ ~ AC (L*,N*) A(L*,N*) B(L*,N») 

pb P,J dX ^ 

where the summations are carried out over the wing planform. For 
the component loadings defined over + he entire planform by 
analytical expressions, the loading is scaled by the appropriate 
power of either semispan or average wing chord, so that the 
resultant component loading lift coefficient is approximately one. 

Interference drag coefficients associated with the configuration- 
dependent loadings are special cases. The configuration- 
dependent loadings are those produced by fuselage upwash, fuselage 
volume asymmetry, and nacelles, as described on page 56. They are 
calculated at discrete points over the wing planform, and 


°D,iJ = C D,J1 ‘IS £ 4C p, i (L * < L * •"*) A(L *’ N,) B(L * 
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TABLE I 

DESCRIPTION OF WING LOADING TERMS 


Loading 

Number Definition 


1. 

2 . 

3. 

4. 

5. 

6 . 

7. 

8 . 
9. 

10 . 

11 . 

12 . 

13 . 

14. 

15 . 

16. 
17 . 


Uniform 

Proportional to x, the distance from the leading edge 
Proportional to y, the distance from the wing centerline 
Proportional to y 
Proportional to x 

Proportional to x(c - x), where c is local chord 

Proportional to x (1.5 c -x) 

Proportional to 2 (1 + 15 A) -0 - 5 
c 

Proportional to (c -x)®® 

Elliptical spanwise, proportional to Vo -y/|) 

Proportional to x, the distance from the leading edge of an 
arbitrarily defined region 

A camber-induced loading proportional to the body bouyancy 
loading 

A camber-induced loading proportional to the body upwash loading 
A camber-induced loading proportional to the nacelle buoyancy 
loading 

The body bouyancy loading 
The body upwash loading 
The nacelle buoyancy loading 



interpolated linearly for a complete definition where needed . 
Since there is no camber surface associated with the fuselage 
loadings, the interference drag coefficient between the "regular" 
loadings and the configuration dependent loadings has only one 
term. To ease the computational task of evaluating component and 
interference characteristics, the configuration-dependent loadings 
are assigned a corresponding wing slope definition of zero. 


The nacelle configuration-dependent loading has yet another 
difference, since the flow fields of the "regular" loadings act on 
the fixed geometry of the nacelles, producing an appreciable axial 
force. The lifting pressure acting on the wing lower surface is 
assumed to propagate along Mach lines from the wing to the nacelle 
area distribution, and the associated buoyancy force is computed. 


O p t im ization of loading. combination. - Given the aerodynamic 

characteristics of the component loadings and their interference 
drag terms, the aerodynamic characteristics of a wing supporting 
any combination of these component loadings can be calculated in 
terms of load strength factors Aj assigned to each. The total 
lift coefficient for n loadings is 


i=n 


C L = 


i=l 


C L,i A i 


(29) 


where c l,i denotes the loading and Aj its load , strength factor. 

The total pitching moment coefficient at zero lift is c 

r ^ mo 


i=n 


C 

mo . 


= V* c A. 

/ r mo. x 


i=l 


(30) 


where 


dC M 

C = r _ £* q 

mo ,i M,i dC T L,i 

Ju 


and the total drag coefficient is 

i=n j=n 


c d = i 


JL - 2 c D,i,j A i A j + 2 c 


i=l 3-1 * 


i=n 

E 

i=l 


DWN ,i A i 


(31) 


c . 

where the terms dwn , i are the axial nacelle force coefficients 

for the component loadings. The contribution of each loading to 

the lifting pressure coefficient at point j can be summed to give 
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(32) 


'P,J 


i= n 


■ £ 


i=l 


'P.i.J 


The load strength factors of the configuration-dependent loadings 
are 1-0, since those loadings are not variable. 


Hence, from the component loadings data, the drag, lift, and 
pitching moment characteristics of an infinite number of wing 
designs may be immediately calculated in terms of the A^ factors. 
The formal optimization of the wing loading involves the 
specification of the set of A^i s which gives least drag subject to 
the imposed constraints. 


According to the Lagrangian method, the solution for the optimum 
A^ values involves the system of linear eguations shown in figure 
4.4-4, where matrix notation is used for compactness. The term 
* c P,i i- s the allowable lifting pressure coefficient corresponding 
to the i th pressure limiting constraint on the wing upper surface 
pressure. 


Partitions in the sguare matrix of the solution have been 
identified by number and letter in figure 4.4-4. The left, 
uppermost partition of this matrix (zone Al) and the adjacent row 
and column correspond to Grant’s original design solution 
(reference 14). The next m rows (zones C1-C5) and columns (zones 
A3-E3) are introduced by m local constraints on the lifting 
pressure coefficient, which are added (if necessary) in applying 
the pressure limiting constraint, as described on page 49. 


The next three rows (zones D1-D5) and columns (A4-E4) correspond 
to the constraints which set the loading factors A. to 1.0 for the 
body bouyancy, body upwash, and nacelle buoyancy loadings. Each 
of these three rows and columns is present only if the 
corresponding configuration dependent loading is used. The right- 
hand column (zones A5-E5) and bottom row (zones E1-E5) correspond 
to the constraint on pitching moment coefficient at zero lift; 
this constraint is also optional. 


The unknowns of the design solution are the n values of the 
loading factors Ai and the (m+5) Lagrange multipliers ( ) . 

These are contained in the left-hand-side column matrix. 
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FIGURE 4.4-4.— OPT! Ml ZA TION MA TRIX 


The five partitions of the right-hand-side column matrix contain: 
either the negative of the nacelle interference drag produced by 
winq lift, or zero if nacelles are not used; the design lift 
coefficient; the constraint values applied to the lifting 
pressure coefficient (these are generated automatically by the 
design computer program, if they are necessary) ; the values of A] 
for the three configuration-dependent loadings; and the design 
pitching moment coefficient at zero lift. 


The maximum number of equations in figure 4.4-4 is 34. This 
number of equations can be reached if 17 loadings are combined, a 
pitching moment constraint is used, and the maximum number of 
constraints on lifting pressure coefficient are required. 


Des ig n camber_sur f ace. - Given the set of Aj factors, the design 

pressure distribution is known and the resulting camber surface is 
calculated using equation 27 on page 40. 


Des ign point opt ions . - Due to options in the wing design program 
to limit the solution extent (e.g., no Cl mo constraint or no 
pressure limiting constraint) , some parts of the solution matrix 
on page 48 may not be used. However, for purposes of describing 
the solution options, it is assumed that all options are chosen. 
In that case, the wing loading solution proceeds as follows, using 
repeated applications of the corresponding matrix: 

(1) The wing loading having minimum lift-dependent drag 
with a constraint on design lift coefficient only, 
is defined. 

(2) A family of wing loadings having minimum lift- 
dependent drag with a constraint on design lift 
coefficient and a series of constraints on 
pitching moment coefficient at zero lift are 
defined. This series of solutions is presented 
in the form of drag due to lift factor, Kjg, 
versus C mo , and is referred to as the "bucket" 
plot. 

(3) The wing loading having minimum lift-dependent 
drag with constraints on design lift coefficient 
and wing upper surface pressure is found. This 
solution may require a cyclic process for the 
pressure constraints that begins with (1) above. 

If (1) satisfies the pressure limits everywhere, 
then this loading is set equal to (1) . If (1) 
does not satisfy the pressure limits everywhere, 
then a pressure constraint is applied at the wing 
planform location where the pressure limit is 
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most strongly violated, and a new solution for 
wing loading is found. Its pressure distribution 
is tested against the limit, and another constraint 
is imposed if the limit is violated. Each cycle 
adds another pressure constraint to the solution 
in addition to the lift coefficient constraint. 

Cycling stops either when a solution load 
distribution everywhere satisfies the pressure 
limit or when the number of solution constraints 
is equal to the number of loadings. 

(4) The wing loading having minimum lift- depen dent 
drag with constraints on design lift coefficient 

and zero-lift pitching moment coefficient is defined. 

(5) The wing loading having minimum lift-dependent drag 
with constraints on design lift coefficient, pitching 
moment coefficient at zero lift, and wing upper 
surface pressure is found. The latter type of 
constraints are imposed, if necessary, in the cyclic 
fashion of (3) . 

Examples of loading solutions (1), (2), (4), and {5) are shown in 

figure 4.4-5. For this case, loading solution’ (1) satisfies the 
wing upper surface pressure constraint, and is therefore identical 
to loading solution {3) . 

The effect of the number of loadings on the "bucket" plot and on 
loading solution (5) is illustrated in figure 4.4-6. Increasing 
the number of loadings from 3 to 7 results in a substantial 
theoretical drag decrease, much greater than the reduction 
obtained in further increase of loadings to 10. 

A discussion of the loading selection and the resultant camber 
surface optimization is presented in part 2 (user's manual). 


Lift Analysis Case 


hoading_f or_a_giyen_camber_su rf ace. - Since the S function for the 
field point element is zero (the element generates no downwash or 
upwash upon itself) , equation 27 can be rearranged to solve for 
the lifting pressure coefficient in terms of the field point slope 
and upstream influences: 


. dz 

AC p (L*,N*) =~g-^ (L*,N*) 


+ — V' R(L* - L,N* -N) A(L,N) B(L,N) C(L,N) AC (L,N) 

'TT £ 4 p 
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Therefore, the lifting pressure distribution, ACp, can be 

determined for a wing of arbitrary surface shape provided the 
calculations are performed in the proper sequence. The order of 

calculation of Acp(L*,N*) is from apex aft, and from the 

centerline to the right hand wing tip for each L* row. In this 

fashion, all pressure, coefficients within the Mach forecone from 
any element will have been previously obtained and no unknown 
pressure coefficients arise in the summation. The element 
fractions A, B, and C are as defined for equation 27. 


Theoretically, ACp (L*, N*) defined by equation 33 is the pressure 
coefficient at the aft mid-point of the L*, N* element. The 
average pressure coefficient for the element, needed in subsequent 
calculations, is calculated by one of two alternative methods: 

1) An approximate average pressure coefficient is 

interpolated from the aft mid-point ACp value and the 
ACp value of the element immediately ahead. (No 
interpolation is performed if L*, N* is a leading edge 
element) . This is the method used in the original 
formulation of the computer program. However, it was 
found that osillations in the ACp values occurred near 
the wing leading edge. So, a 9 point smoothing equation 
was applied after unsmoothed pressures were computed for 
the entire wing, which essentially removed the 

oscillations and brought the smoothed pressure 

distribution into good agreement with established 
results from other theoretical methods. The smoothing 
equation is of the form 

(34) 

0.2A(L*-4)AC (L*-U) 0.1*A(L*-3)AC (L»-3) + 0.6a(L«-2)AC (L»-2) 

P P P 

♦ 0.8A(L*-l)AC (L*-l) + A(L*)AC (L») + 0.8A(L* +l)AC p (L* + l) 

+ 0.6A(L»+2)AC p (L«+2) + 0.1»A(L»+3)AC p (L*+3) + 0.2A(L*+U)AC p (L*+U) 

AC p (L # ) " 0.2A(L*-U) + 0.1*A(L»-3) + 0.6 a(L»- 2) + 0.8A(L»-l) + A(L») ♦ 0.8A(L»+l) 

+ 0.6A(L*+2) + O.UA(L*+3) + 0.2A(L»+4) 

where the A values account for element fractions. A is 
zero if the corresponding element is ahead of the wing 
leading edge, or aft of trailing edge in the case of a 
subsonic trailing edge. The supersonic trailing edge 
solution is treated as a special case, and the trailing 
edge is extended four elements to provide pressure 
coefficient data to fill out the smoothing equation. 

2) The other average ACp method is an improved technique, 
which uses an aft element sensing approach. This method 
involves solving for preliminary ACp results for a given 
L*, N* element and the element immediately aft, then 
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following it up with a second calculation to refine the 
preliminary results. The procedure is detailed in 
reference 12, but is summarized below and in figure 4.4- 
7. 


a) 


b) 


c) 


Calculate preliminary Acp values for a given L* = 
constant row. Designate as Ac a (L*, N*) 


Calculate preliminary Acp 
immediately aft = L’i' + I 

Designate as AC , (L# N*). 

P.b 


values for the element row 

, using AC values. 
p,a 


Finalize Cp values for the original L* = constant 
row from one of the two following equations. 

For leading edge elements: 


AC (L*,N*) ±- 


_ 1 _ } -j . A(L* ,N* ) 

2 [ 1 + A(L* ,N* ) 

1_ \ A(L* ,N*) 

2 [l + A(L*,N») 


AC (L*,N«) 
P » a 


( 35 ) 


AC ,(L*,N«) 
p.b 


For all other elements: 


AC 

P 



(L*,N*) + -r- AC . (L*,N*) 
4 p.b 


( 36 ) 


The aft element sensing technique produces a substantial reduction 
in pressure coefficient oscillation. In cambered wing test cases 
run with the method, however, some oscillation was found to 
persist. Therefore, a 3 term smoothing equation was added (after 
all wing pressures are calculated) , of the form: 


AC (L*)= 
P 


0.5A(L*-l)AC (L*-l) + AC (L* ) + 0.5A(L*+l)AC (L*+l) 

P P P 

0.5A(L»-1) + 1.0 + 0 . 5A ( L*+l ) 


( 37 ) 


In the computer program one of the two alternate smoothing 
methods is selected by an input code. Both produce essentially 
the same answers . When the pressure limiting option is used 
(discussed later) , the aft-sensing smoothing technique is 
automatically selected. 


For the analysis solution, the pressure coefficients for all 
elements must be calculated. The force coefficients are, 
therefore, calculated from direct summations of local pressures 
applied to each element, rather than employing a spanwise 
integration as in the design case. Lengthwise and spanwise lift 
distributions are obtained by summing the lift in the 
corresponding element rows. 
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FIGURE 4.4-7.-APPL/CAT/0N OF AFT-ELEMENT SMOOTHING TECHNIQUE ( FROM REF, 12) 



No "leading-edge-suction" force is included in the drag summation, 
which accounts only for the lifting pressure coefficients acting 
on the local wing slopes. 

Analysis of complete co nfi g ur ati on . - preceding discussion 

describes the theory used in the lifting pressure calculations for 
a wing camber surface at a selected angle of attack. The analysis 
program may also be used to calculate lifting pressure 
distributions and force coefficients for complete configurations 
over a range of angles of attack, adding in the effects of 
fuselage, nacelle, canard and/or horizontal tail, as applicable. 


The program actually carries two solutions along: one for the 
configuration at its input angle of attack, the other the 
incremental solution per degree angle of attack (called the flat 
wing solution). The interference terms associated with the two 
solutions acting on the other surface shape are also calculated. 
The summation of these effects into the drag polar and other force 
coefficients is performed by superposition, as described on page 
64. 


Calculation of the complete configuration lifting pressure 
solution involves up to 7 principal tasks: 

• Isolated fuselage upwash field 

• Nacelle pressure field acting on wing 

• Pressure field due to asymmetrical fuselage volume 

• Wing/canard solution in presence of fuselage upwash 
field. Effects of the wing pressure field acting 
on nacelles, are also calculated 

• Fuselage lift distribution in presence of wing 
downwash field 

• Horizontal tail solution in presence of fuselage 
and wing flow fields. 

• Superposition of solutions 

The calculation logic of these tasks is described in the following 
sections: 


Isolated^. fuselage upwas h fi eld. - The fuselage upwash (or 

downwash) field in the plane of the wing, canard or horizontal 
tail is obtained from slender body theory. From reference 15, the 
velocity potential of a slender body of revolution, shown in 
figure 4.4-8, is: 

,x 


(38) 


1 A( Qi _ 1 dA(x) „„ 8r 13 

TT* U * r ’ e) - Tit ST- ln ~r -■&-& 


r 


dA(x) 

“dx 


£n(x-x 1 )dx 1 + [cx+e(x)]R (x) 


cos 6 


The nondimensional radial velocity component is 
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Mach "i 
line 


Upwash propagates along Mach lines. 
Upwash at arbitrary point P is: 


w o = u = JL dA(xj z_ _ { } R 2 (x) /z 2 -y 2 \ 

~ ° 2vr dx r 2 r 2 r 2 ) 

W. = u - - R 2 (x) /z 2 -v 2 \ 

U r 2 V r 2 / 

Where: 

u Q = upwash at input fuselage incidence (a= 0) 

u = upwash increment per unit a 

A(x), R(x) = Fuselage cross-sectional area, radius 

e (x) = Fuselage incidence due to camber 
(relative to reference plane) 

z,y,r = location of point at which 
upwash is calculated 

a = angle of attack of reference plane 
U = free stream velocity 

FIGURE 4.4-8. —FUSE LAGE UPWASH CALCULATION 
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JU = — L ^AU) _ c L +_|00 r 2 (x) cos e (39) 

U ? r 27rr dx r 2 

and the tangential component is 

2 

-jj <t> e = -[a + e(x) J ■ " g X ) Sin 8 (1+0) 

r 

The vertical velocity component (upwash angle) is 

i = i* r co S e - i * e sin e 0*1) 


Substituting for 4> r and <f> 0 and converting 
gives 


w _ JL 
U 2 it 


dA(x) 

dx 


- [a + e(x) ] 


R 2 (x) 

2 

r 


to cartesian 



coordinates 

(42) 


When the fuselage is at input incidence, the fuselage upwash is 

I J ■ I (1*3) 


u = i- “hd ’ . e(x) “fw 

r r \ r 


o 2ir dx 

and the incremental upwash per degree angle of attack increase is 


• d w/u R 2 (x) / z 2 - y 2 
U= da r 2 r 2 


(44) 


Slender body theory assumes that the upwash field propagates 
normal to the body centerline. In the analysis program, the 
fuselage upwash is propagated along Mach lines, as reguired by the 
characteristics of supersonic flow. 


Na celle pressure_field_a.ct.ing on wing . - The pressure field caused 
by the nacelles on the wing is calculated by the Whitham solution, 
as described in Section 4.3. The nacelles are assumed to be 
bodies of revolution, and the pressure signature due to each 
nacelle is calculated for a series of spanwise stations from wing 
root to tip. The composite pressure signature at a given spanwise 
station is the sum of the individual nacelle signatures, with all 
pressure coefficients doubled to account for reflection from the 
wing. Either ’’wrap” or ••glance” nacelles may be used, as 
discussed in Section 4.3. 


Pressure field due t o asymmetric fuselage v olu me. - Another term 

in the superposition approach is associated with the growth in 
fuselage area above and below the wing. If the area growth is 
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asymmetric . (e.g., a low wing configuration), a differential 
pressure across the wing plane is created which is additive to the 
effect created by the fuselage upvash field. 


an approximate method of computing the asymmetric fuselage 
pressure field is used, using the area distributions shown in 
figure 4.4-9. To the actual forebody area distribution is added 
the growth in fuselage area in the wing region relative to the 
area at the leading edge of the wing-fuselage intersection. This 
is done for both the below- amd above-wing area distributions. 
Pressure signatures for both body representations are computed 
using the Whitham technique at a series of spanwise stations, 
doubled to account for reflection from the wing, and differenced 
to get the lifting pressure distribution due to asymmetric 
fuselage volume. 


If the wing trailing edge is subsonic, the pressure fields above 
or below the wing revert to that corresponding to the actual 
fuselage area distribution in the region aft of a Mach line from 
the trailing edge of the wing-fuselage intersection. 


Use of the asymmetric fuselage calculation is controlled by an 
input code in the computer program. If the asymmetric calculation 
is not requested, the fuselage pressure field corresponding to a 
mid-wing arrangement is calculated, so that there will always be 
a thickness pressure field due to the fuselage, if present, for 
use in limiting pressure calculations (discussed later) . 


Wing/Canard Solution in Presence of Fuselage and Nacelles 


The wing or canard lifting pressure solution in the presence of 
fuselage and/or nacelles is performed as described for the 
isolated wing, with the following alterations: 

• The local surface angles of attack are increased by the 
fuselage upwash values for the purposes of computing 
pressure coefficients. (To compute drag, the pressure 
coefficient is applied to the local surface slopes, 
only) . 

• In the region of the wing covered by the fuselage, the 
wing slopes. are zeroed. This reflects the fact that the 
wing lift distribution in this region is of the "carry- 
over" type, only. (To compute drag, the carry-over lift 
is applied to the fuselage camber line shape) 
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FIGURE 4.4-9.- FUSELAGE AREA REPRESENTATION FOR ASYMMETRIC 

FUSELAGE BUOYANCY FIELD 
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• The canard pressure distributions are computed in the 
same fashion as the wing pressure distribution, with the 
canard also creating a downwash (or upwash) field for 
the wing. 


The nacelle or asymmetric fuselage volume pressure fields are 
superimposed upon the wing elements to obtain the effect of the 
nacelles or fuselage asymmetry on the wing. If nacelles are 
present, the buoyancy drag of the wing lifting pressures acting on 
the nacelle area distribution is also calculated. This term is 
obtained by transferring local wing pressures aft along Mach lines 
to act upon the nacelle area growth, in a fashion paralleling the 
wing- thickness-on-nacelle term of the near-field wave drag 
program. 


Fusela ge lift distribution in presence of win g downw as h fi e ld . - 
From slender body theory , the fuselage lift distribution is given 
by 

c g • R(x) = [a + e(x)] + -|^ [e(x) • A(x)] (1*5) 


using the notation of figure 4.4-8. 

The first term of the equation is the slender body term associated 
with a straight body at angle of attack; the second term is due to 
the curvature change associated with body camber. 


In computing the fuselage lift distribution, the downwash (or 
upwash) from the wing/canard is included in the local angle of 
attack term, e . For drag computations, the lift is applied only 
to the physical angle of attack value, a + e{x). 


Fuselage force coefficients (lift, drag, and pitching moment) are 
calculated from the slender body lift distribution and converted 
to a wing reference area basis for summation into the complete 
configuration characteristics. 


Horizontal tail so lu tion. - The horizontal tail lifting pressure 
distribution is calculated in the presence of the fuselage upwash 
field and wing/canard downwash field. The solution employs the 
same logic of partitioning the tail into exposed and carry-over 
regions used in the wing and canard pressure computations. 
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Since the horizontal tail is assumed to be aft of the wing, the 
fuselage and wing/canard upwash (or downwash) field is first 
calculated, then the tail lift distribution and force coefficients 
are computed for the desired tail angles of attack. These are 
summed into the wing/fuselage, etc., solution for each tail angle, 
resulting in a set of force coefficients for each tail setting in 
the same fashion as wind tunnel data are obtained. 


Calculation of the wing/canard downwash field is performed by 
extending the wing grid system aft to include the horizontal tail. 
With all pressure coefficients on the wing or canard previously 
established, the effect of the wing/canard on the horizontal tail 
is obtained from 

Wing __ 

A w c = Z R (L» - L, N* - N) A (L, N) B(L, N) C(L, N) AC (L, N) 

* I ' (46) 

Canard __ 

+ 53 R (L* - L, N*- N) A(L, N) B(L, N) C(L, N) AC (L, N) 

T ^ 

in which the wing, canard, and tail are all assumed to be located 
in essentially the same plane. 


Calculation of the horizontal tail lifting pressures are then 
performed using the analysis form of the lifting pressure eguation 


AC p (L*, N*) = -jj (L* , N*) 

(kl) 

. Tail _ 

+ - 23 R (L* - L, N* - N) A(L , N) B(L, N) C(L, N) AC (L, N) + A u _ 
it tr p w,c 


with the wing/canard term added to the summation of upstream 
effects. As in the wing/canard pressure solution, the fuselage 
upwash is added to the physical tail slope for the purposes of 
computing lifting pressures, but removed for the calculation of 
the drag term. The carry-over lift of the horizontal tail is 
applied to the fuselage camber line to obtain the drag 
interference of the tail on the fuselage. 


The effect of wing downwash on the horizontal tail lifting 
pressures is significant. The data of figure 4.4-10 show the 
theoretical drag and pitching moment increments at constant total 
lift for a typical configuration, with and without wing downwash 
included. (Comparisons of the theoretical calculations with 
experimental data are presented in the user's manual, part 2.) 
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Superposition of Solutions 


The arrangement of terms in the superposition of the complete 
configuration force coefficients is: 


°l =£ c l flic, 

a = 0 


C = Sc 
m m 


a = 0 


+ a S 


m 


(48) 

(49) 


C D * 2 °D . + K 1 (C I - C L . > + K 2 (C L - C L . >' 

a = 0 • a = 0 a - 0 


(50) 


where the lift, drag, and pitching moment coefficients at «= 0 
correspond to the individual fuselage, wing, tail, etc. , solutions 
at input configuration incidence, and the AC L and Ac^ 
coefficients correspond to the incremental (flat wing) solutions 
per degree angle of attack. 


The drag coefficient equation includes the interference terms 
between the solution at input incidence and the incremental 
solution per degree. With subscripts C and F denoting input 
incidence and the flat wing terms, respectively, the drag 
coefficient equation is: 


c ° = % * %c 


+ (L ) — + -2- C_ a 
d cof “f “f l f 


(C T 


+ C T 


) 


0.01745 (C T - 


or CL = 


V 


FOC 


COF 


- y ♦ 


C L >' 

L c 


(51) 

(52) 


where a= configuration angle of attack 


a = flat wing incidence. 



.01745 radian 




a e “c 





with a e = 



local element area corresponding to AC 
reference area 
local element slope 
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interference drag of flat wing 
pressures acting on cambered 
wing slopes, per degree. 





interference drag of cambered 
wing pressures acting on flat 
wing incidence, per degree. 


with the summations carried out over the wing planform. 


Therefore, for each of the solutions involved in the 
superposition, it is necessary to calculate both an input 
incidence solution, a flat wing solution, and the interference 
drag terms between the two solutions. All of the force 
coefficients are referred to an input area, moment center, and 
moment reference length in the summations. 


Limj ti r.q pressure solution. - Linear theory imposes no 

restrictions on the allowable wing pressure coefficients. 
Particularly in the case of subsonic leading edge wings, large 
upper surface pressure coefficients may be computed near the 
outboard wing leading edge that exceed (are more negative than) 
pressure coefficients corresponding to a vacuum. These pressure 
coefficients can occur at moderate wing angles of attack, due to 
the strength of the upwash generated by the inboard wing. 


Several investigators have examined experimental data to define 
the minimum attainable upper surface pressure coefficients. One 
such correlation is shown in figure 4.4-11 (from reference 16), 
which indicates a physical pressure coefficient limit of 
approximately .8 vacuum. The mechanism of this limiting is 
associated with local leading edge flow separation. However, the 
phenomenon may be approximated in the wing analysis program by 
limiting the calculated linear theory pressure coefficients to a 
specified fraction of vacuum C p . 


In the analysis program, pressure limiting is an optional feature, 
controlled by an input code. If limiting is requested, then a set 
of configuration angles of attack is required for the solution, 
since superposition will not apply after limiting occurs. Also, 
a definition of the wing thickness pressures is required 
(transferred over from the near field wave drag program module) , 
since it is the sum of lifting plus thickness pressures that is 
limited, rather than the lifting pressure alone. 


Pressure coefficients, as calculated, are separated into upper and 
lower surface values. Thickness pressures (wing plus fuselage) 
are added and the upper surface pressure coefficient tested 
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FIGURE 4.4-11. -MINIMUM ATTAINABLE PRESSURE COEFFICIENTS (FROM REF. 16) 




against the limit pressure coefficient. If the limit is exceeded, 
the calculated lifting pressure coefficient is adjusted to a level 
which causes the total upper surface pressure coefficient to match 
the limit value. 


Comparisons of the "limiter” calculations with measured force and 
pressure coefficient data are presented in figures 4.4-12 through 
4.4-16 for 2 arrow wings of references 1 and 2. Both were 70° 
sweptback wings, having 3 percent biconvex airfoils, and wing 
design lift coefficients of 0 (flat) and .08. 


The basic thickness pressure comparison is shown in figure 4.4-12. 
Comparisons of upper and lower surface pressure coefficients at 
lifting conditions of C L « .13 and .25 are shown in figures 4.4-13 
and 4.4-14. In the case of the flat wing (figure 3.4.12), several 
fractions of vacuum pressure coefficient are illustrated; the C L 
= .08 wing data are shown for .7 vacuum limiting only. 


The limiter feature exhibits considerable improvement in the 
detail pressure coefficient comparisons over unlimited linear 
theory as C L increases. 

Force coefficient comparisons for the two wings are shown in 
figures 4.4-15 and 4.4-16. Pressure limiting at approximately .7 
vacuum improves the comparisons appreciably at the higher C L 
values . 
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FIGURE 4.4-13.-PRESSURE COEFFICIENT COMPARISON, WING 1 (FLAT WING) 





I 


FIGURE 4.4-14.-PRESSURE COEFFICIENT COMPARISON, WING 2 

(TWISTED AND CAMBERED WING, DESIGN C, = .08) 











FIGURE 4.4- 16 - TEST- THEOR Y COMP A RISON, WING 2 (C L 








5.0 INTERACTIVE GRAPHICS 


Interactive graphics for use with the design and analysis system 
are optional, employing the NASA- LRC cathode ray tube (CRT) 
display and associated software. Use of the graphics option is 
requested by an executive control card (described in part 2, 
user' s manual) . 

The principal uses of the graphics routines in the design and 
analysis system are to display the configuration, edit input 
geometry, and to display and/or alter the basic program 
calculations. Limited capability to redirect the system 
calculation sequence is available from the CRT console during 
program execution. 

The CRT display and program coding for the interactive graphics 
setup are based on the NASA-LRC system. However, all display 
portions of the coding are subroutined or overlaid from the basic 
programs, so that the system could be readily converted to other 
CRT arrangements. 

Details of the interactive v graphics portion of the design and 
analysis system are given in part 2 (user's manual). 
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